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We simulate electron transport through graphene nanoribbons of experimentally realizable size 
(length L up to 2/im, width W ^ 40nm) in the presence of scattering at rough edges. Our numer- 
ical approach is based on a modular recursive Green's function technique that features sub-linear 
scaling of the computational effort with L. We investigate backscattering at edge defects: Fourier 
spectroscopy of individual scattering states allows us to disentangle inter- valley and intra- valley scat- 
tering. We observe Anderson localization with a well-defined exponential decay over 10 orders of 
magnitude in amplitude. We determine the corresponding localization length for different strength 
and shape of edge roughness. 
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FIG. 1. (a) Graphene hexagonal lattice with lattice constant a — 1.4A. The unit cell (shaded area) contains two carbon 
atoms A and B belonging to the two triangular sublattices connected to each other by the displacement vector tab [see inset] . 
Each atom has three nearest neighbours [smallest (blue) circle], six next-nearest neighbours [medium (green) circle] and three 
second-nearest neighbours [largest (red) circle] . (b) The conical dispersion relation with trigonal warping of an infinite graphene 
plane near the K and K' points, as obtained by the third-nearest neighbour tight-binding approach. 



I. INTRODUCTION 

The experimental realization of graphene, i.e., of a monolayer of carbon atom^^H^has opened up a rapidly developing 
field of fundamental and applied physics. The topology of theplanar honeycomb lattice [figure [TJa)] with the resulting 
peculiar band structure near the K and K' points [figure [T'b), for a review, see"^ ^] gives rise to many novel and 
intriguing physical properties, including the room temperature quantum Hall effect, minimum conductivity at the 
Fermi energy as well as possible applications for spintronics. Recent advances in fabricating width- modulated graphene 
nanoribbons helped to overcome intrinsic difficulties in creating tunnelling barriers and confining electrons in graphene, 
where transport is dominated by Klein tunneling-related phenomenaPEI. Graphene quantum dots have been fabricated 
and Coulomb blockade^ quantum confinement^^ and charge detectiorP^ have been demonstrated. 

The electronic properties of the perfect honeycomb lattice are meanwhile theoretically well understood^. However, 
in realistic graphene devices finite-size effects and imperfections play an essential role, especially for transport through 
confined structures. The importance of such effects results from the gapless band structure of graphene which does 
not allow straightforward confinement by electrostatic potentials. Devices thus have to be cut or etched resulting in 
rough edges. In turn, properties of the ideal graphene band structure cannot be invoked when simulating quantum 
transport through realistic devices in the presence of randomly shaped boundaries. Moreover, recent results underline 
that an incoherent Boltzmann transport approach, unlike a full quantum- mechanical calculation, fails to reproduce 
experimentally observed conductance signatures for impurity scattering^^. However, the application of numerical 
methods for a full quantum mechanical simulation of graphene ribbons of realistic size constitutes a considerable 
challenge. A method of choice is the widely-used recursive Green's function technique^^ which is well suited to 
treat scattering structures such as wires or ribbons extended in one of their dimensions and usually implemented 
within the Landauer-Biittiker frameworlP^ for calculating transport coefficients. This technique has meanwhile been 
incorporated successfully into both an equilibriunF^*^ and a non-equilibr iuni^ (Keldysh) description. Several 
variants of this method have been put forward employing specific symmetries of the system^^ if applicable, or using 
recursive algorithm^^^M^. In this article we present an extension of the modular recursive Green's function method 
(MRGM)'^ designed to treat graphene nano-ribbons with edge or bulk disorder. We address disorder scattering in 
graphene nanoribbons both on the microscopic level of specific lattice defects as well as on the macroscopic level 
of current measurements for which, as we will demonstrate, the details of the underlying disorder scattering play a 
crucial role. 

This paper is organized as follows: We briefly review key properties of the band structure of "ideal" infinitely 
extended graphene in the absence of disorder in section |TI[ In section [TTTj we introduce the application of the MRGM 
to finite-size graphene structures which allows us to treat extended structures efficiently due to the favorable scaling 
of the numerical effort with the linear dimensions of the ribbon. Applications to transport through rough-edged 



graphene nanoribbons will be presented in section IV followed by a short summary (section [V] 



II. TIGHT-BINDING SIMULATION OF GRAPHENE BAND STRUCTURE 



The ideal, infinitely extended graphene sheet features a honeycomb lattice made up of two (A and B) interleaved 
triangular sublattices. It can be described in tight-binding (TB) approximation by the HamiltoniarP^ 

H = ^ Vi {(l)i^s \ - ^ 7i,j + h.c. , (1) 
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where the sum extends over pairs of lattice sites, |^j,s) is the tight-binding orbital with spin s at lattice site 
j, Vi is a locally varying potential (onsite energy), and jij is the hopping matrix element between lattice sites i 
and j. Within our TB approximation, we include third- nearest-neighbour coupling [see figure [TJ a)] using orthogonal 
tight-binding orbitals. This allows for four free parameters, namely the site-energy eo and the overlap integrals 7^, 
i = 1,2,3, representing the interaction with the first, second and third nearest neighbour, respectively. We choose 
the 7i by fitting to ab-initio calculations, taken from Reich et aLl^^M^, arriving at 71 = —3.145, 72 = —0.042, and 
73 = -0.35. 

The dispersion near the non-equivalent K and points [figure [TJb)] resulting from the diagonalization of equation 
([1]) features [for large distances k from K{K')] deviations from a perfect cone refiecting the infiuence of the hexagonal 
lattice. The cone becomes squeezed along the K — K' directions, an effect known as trigonal warping^ Near the 
K point and for small k the band structure of equation ([T]) can be approximated (assuming that 1^ <C t) by a conical 
dispersion relation around the K poinl!^, 

E{k + kK) = E{kK) + kdkE{kK) + 0{kj^) ^ VFh\kl (2) 

where we have set E{kK) = 0. Note that the above expansion ignores both the length scale of the graphene lattice 
constant a = 1.4 A and the preferred directions of the lattice due to the discrete lattice symmetry. In this low- A: limit, 
the Hamiltonian, equation ([T]), can be approximated by the Dirac Hamiltonian, 

H = Vf {pxCTx (8) To + PyCFy (g) Tz) (3) 

with a and r being the Pauli spin matrices acting on the pseudo-spin and valley degrees of freedom. Analytic solutions 
for an infinitely extended graphene sheet described by equation (|3| yields plane waves \k) where the angle of the k 
vector 6>/c, 

Ok = tan-^{ky/k^), (4) 

connects relative amplitudes on the A and B sublatticd^l, 

|k) = e^^-'^(|A)+e^^^ \B))/V2. (5) 

Consequently, the pseudo-spin projection along the direction of propagation, the "helicity" , 

h=ia ■ k)/ |fc| (6) 

is conserved refiecting the chiral symmetry of the ideal graphene sheet in the \ow-k limit. The additional degeneracy of 
two non-equivalent cones ( "valleys" ) at the K and points in the reciprocal lattice allows to formally represent the 
low-energy band structure near E = in terms of Dirac-like four-spinors {ip) = {ip^ , , ip^', ip^') with amplitudes 
for the A — B sublattice in real space and K — K' in reciprocal space. The sign of 9^ [equation Q] is reversed upon 
transition from K io K' . (Note that physical spin is not included in the present analysis.) 

This Dirac-like picture may serve as valuable starting point for the analysis of finite-size and edge effects on transport 
incorporated within the tight-binding Hamiltonian [equation ([T])]. Chirality [equation (|6|] is preserved in the presence 
of slowly (on the scale of the C-C bondlength) varying perturbations, suppressing backscattering^ 

P(k ^ k') = Kk'l V |k)|' oc cos2[(0k - ^k')/2]. (7) 

Conversely, rough edges may break chirality resulting in non- vanishing backscattering on the same cone and, at the 
same time, in coupling of the K and K' cones. In the following, we will investigate in detail the effects of short-range 
defects^^ (e.g. edges) and deviations from the continuum Dirac picture on the transport properties of rough-edged 
graphene nanoribbons. 

III. NUMERICAL METHOD 

For the numerical treatment of finite-size graphene fiakes and ribbons a number of simulation algorithms have 
meanwhile been progosedSElEil the following an extension of the Modular Recursive Green's function 

Method (MRGM]PE33I1 applied to the third-order TB Hamiltonian [equation The key idea of the MRGM is to 
break down a large device into independent smaller modules, each of which can be computed efficiently [see figure [2]. 
The Green's functions Gu of the different modules with width W and length L are then combined to the desired device 
geometry using a small number of Dyson equations. In this Article, we introduce an efficient method to calculate the 



4 



(a) . (b) 



I 



(d) 



I : I I I I I I 
I I I 1 I I I 1 
I : I I I I I I 

V.V.ViV. 




:i I I I I I 

: I I I I I 

:i I I I I I :l I I I I I 

:.V,V.".V.ViV. 

: 1 I I I I I: I I I I I 
■I I I : I I ■! I I I I I 

I I I I I !■ 

: I I I I I : 



(c) 




I I I I 1 I 
I I I I I I 
I I I I I I 
I I I I I : 
I I I I I I 

<' ' r ' ' > 



■I I 1 1 I I 1 I : 
I I I I : I I 1 
:i 1 I I 1 1 1 I : 



Wu 





FIG. 2. Assembling a rough-edged nanoribbon by combining several modules: (a) a chain of carbon atoms in transverse direction 
(y) is periodically repeated to yield (b) a half-infinite graphene ribbon along the x direction, (c) Using Dyson equations [see 
equation (15)], a rectangular region can be separated from the half-infinite ribbon, (d) A rough-edged ribbon can now be 
assembled by randomly combining rectangles of variable length Lu and width Wu. (e) For an arbitrarily shaped graphene 
scattering structure, the modular approach can still be used, although the Green's function of each module has to be calculated 
by direct inversion. 



Gu'- the associated numerical effort becomes independent of module length L. The algorithm involves the somewhat 
counter-intuitive steps to first calculate infinitely and semi-infinitely extended ribbons, i.e., modules of width Wu but 
1/ = oo, from which rectangular modules of finite length Lu are "cut out" as needed by applying the Dyson equation 
"in reverse". The obvious advantage of this approach is that the computational effort becomes independent of Lu- 
This approach is particularly advantageous for the simulation of weakly disordered graphene nano-ribbons: for weak 
disorder the spacing between individual defects both in the bulk and at the edges of the ribbon is large as compared 
to the lattice spacing. We can thus simulate the region between neighbouring defects by a single graphene module 
with perfect boundaries and place adjacent defects at the module boundaries [see red dotted lines in figure [2|d)]. This 
efficient calculation of the extended rectangular modules in between two defects is key to simulating large devices 
with length Ltotai up to several micrometers (or ^ 10^ hexagons in one direction). 

As a prototypical example, we build up an infinitely long nanoribbon with ideal zigzag boundaries (along the 
X direction) by periodic repetition of a chain of carbon atoms of width W in transverse {y) direction [see figure 
[2|a)]. Other geometries and boundaries can be treated analogously. The Hamiltonian H of the ribbon can thus be 
decomposed into a matrix Hq describing the Hamiltonian of the vertical chain, and the coupling matrix Hj describing 
the connection between two adjacent chains^^. 



H 



Ho Hj 

Hj Hq Hj 

Hj Ho 



(8) 



The solution of the Schrodinger equation for the infinite ribbon can be written in terms of an ansatz for a Bloch wave 



(9) 



with \xn) the transverse eigenfunction. For expanding the Green's function in terms of Xn we need a complete set 
of transverse eigenfunctions including all evanescent modes in the sum [equation ([9|] and the subsequent equations. 
The resulting generalized eigenvalue problem for e*^^^ and \xn) gives n left (right) -moving states \xj) {\Xj))i with 



corresponding momentum 
\x,)e"'^^{Xj\ iDj{x) = \xj) 



kj 



{kj) in X direction. In the following, we introduce the shorthand notation Dj{x) 



(Xjl) for the projections onto the right (left) moving Bloch states. From the Bloch 
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(10) 



states the Green's function of the infinite ribbon follows as^ 

f T.^=iDj{x-x')V-\ x>x' 
G^(x,x') = I 

[1 

with the hopping matrix 

AT 

V = Y,Hl[Dj{-Ax)-Dj{-Ax)]. (11) 
i=i 

The Green's function of the half-infinite ribbon Gl(Gr) extending from xo to — oo (or +00) can be written as 

Gr,l{x, x') = G^{x, x') + G%^{x, x') (12) 

with 

Gl{x,x') = ^Dj{x - xo)Dj{xo - x')V-\ (13) 

j 

G%{x,x') = ^Dj{x - xo)Dj{xo - x')V-\ (14) 
j 

satisfying the boundary conditions Gr^l{x^x') = for ah x or x' located at the end of the half-infinite ribbon 
{x^x' = xq). For an intuitive interpretation of equation (12) consider a disturbance from a point source at x' . It 
reaches x by two paths: the direct propagation from x' to x, given by the Green's function of the infinite ribbon, and 



the propagation from x to xq [given by D^{x — xq) in equation (14)], where the wave is reflected at the end of the 
ribbon and then propagates from xq to x' [given by Dj{xo — x') in equation (14)]. Note that the numerical effort to 
calculate and Gr,l is controlled by the transverse width of the ribbon Wu and the number of transverse modes 
\Xn) to be included while the x-dependence is given analytically. This scaling behaviour is key to calculate Gu for 
the rectangular ribbon of arbitrary length Lu by solving the Dyson equation 

Gl,r = Gu + GuHiGl.r (15) 

in reverse for Gu instead of for Gl.r- Consequently, the numerical effort to calculate Gu is independent of Lw In 
the final step, a rough-edged nanoribbon can now be assembled by successively joining rectangular ribbons Gu^ of 
varying length Lu and width Wu (average width W = 60nm) using the Dyson equation in forward direction, 

G = G^^^ ^G^^^HiG. (16) 



In our simulations, we treat ribbon lengths of several micrometers, and average over 100 random realizations ^ of 
edge roughness to eliminate non-generic features of particular ribbon configurations. To assemble such very long 
disordered ribbons, we start with a set of Nb different modules Mi, . . . , Mat^ and combine them to obtain a larger 
module Mnb-\-1' We connect the calculated modules in a random permutation V, M^b+i = V{Mi + . . .-\-Mnb) (e.g., 
for Nb = 5, Me = M3 + Ml + M5 + M4 + M2). This procedure is repeated iteratively [i.e., M7 = V{M2 + . . . + Me), 
formally equivalent to the composition rule of a Fibonacci sequence], creating an exponentially growing, pseudo- 
random sequence of modules. The interfaces between modules that include the disorder are randomly determined at 
each iteration step to avoid periodic repetition. 

If a more general shape of the scattering geometry is desired [e.g., for a non-separable disorder potential or curved 
boundaries, see figure [2|e)], the partitioning into modules and the subsequent efficient buildup of long structures is 
still readily possible. Only the first step of our algorithm has to be modified. The Green's function of individual 
modules is directly calculated by inversion, i.e., G = {E — H)~\ using, e.g., a parallelized sparse-matrix solveM 
Subsequent application of Dyson equations allows to assemble complex scattering geometries. 

We find that the computing time r of our approach scales as ~ due to the cubic dependence of the eigenvalue 
problem and of the matrix multiplications, r scales linearly with the number of building blocks Nb used to set up 
the geometry, and logarithmically with L/Lw Numerically, we find for r, 

T[ms\ - a A/'b In (^-^^ • W^ (17) 
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FIG. 3. (a) Eight building blocks of different length / zb A/ = 3 ± 2nm and height in a range W ± AW/ 2 were used to assemble 
(b) a rough-edged graphene nanoribbon (different shades of grey for clarity), (c) Ensemble- averaged conductance G of 40nm 
wide graphene ribbons of length L = lOOnm with different amplitude of edge roughness AW as a function of back-gate voltage 
Vbg- The conductance of a ribbon with perfect zigzag boundaries is shown as dashed black line. Arrows (t) mark dips in the 
conductance (see text). The shaded area highlights the voltage interval of increased conductance in the ideal ribbon due to 
states localized at the zigzag edge (see text). The solid A [open A] triangle marks the back gate voltage of individual scattering 
states displayed in figure |4] (a) [(b)]. (d) Dispersion relation k[E{V)] of an ideal 40 nm wide graphene zigzag ribbon, enlarged 
around the K and K' points. 



with W given in nm. We have determined a prefactor a ~ 5 from calculating the full scattering problem averaged 
over 100 configurations, when computing on 3 AMD Opteron processors (24 cores) at 2.2 GHz. Glearly, the prefactor 
strongly depends on the details of the employed hardware (i.e., network speed, cache size, compilation flags, etc.). 



while the scaling (17) does not. 

We note that the application of the algorithm presented here is not restricted to graphene nanostructures: any 
modular scattering system, which is build from modules along the lines of figure 2 (a)-(c) can be treated analo- 
gously. Possible applications include acoustic cavities, conventional semiconductors, topological insulators, or neutron 
scattering devices. 



IV. RESULTS 



A. Transport coefficients 



For conventional semiconductor heterostructures (e.g., quantum dots made of GaAs-AlGaAs), confinement is usually 
achieved by electrostatic gates resulting in smooth dot boundaries. Such confinement is not realizable for graphene 
due to its gapless band structure. While several theoretical concepts for opening a band gap have been proposed, the 
majority of experiments have achieved confinement by patterning of graphene nanodevices with oxygen plasma etching. 
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FIG. 4. Scattering states of rough-edged graphene nanoribbons at selected back gate voltage: (a)-(c) Vbg = — 5 V [corre- 
sponding to open triangle in figure [sj^c)], (d)-(f) Vbg = —15 V [solid triangle in figure [s] (c)]. Panels (a), (d) show the entire 
scattering wavefunction, while panels (b,e) [(c,f)] feature projections onto the A [B] sublattice respectively. Frames to the right 
show zoom-ins of wave function enhancements at upper (lower) corners marked by red arrows in (d)-(f), the positions of the 
carbon atoms are marked by white dots as guide to the eye. 

chemical vapor deposition, specially prepared SiC substrateP^, or chemical etching. These techniques, however, do 
not result in well-defined armchair or zigzag edges b ut in a rough-edge pattern featuring armchair and zigzag elements 
as well as adsorbates at the dangling carbon bonds^ leading to an irregular edge structure. Edge effects can thus 
be expected to strongly influence the properties of graphene nanodevices. 

We simulate the influence of edge scattering on transport through graphene nanoribbons by randomly varying the 
widths of the rectangular modules which build up the ribbon in the range W = 40±1 nm [see figure [sja)]. Numerical 
tests show that a random sequence of Nb = 5 different module widths represents a good compromise between a high 
degree of randomness and limited computational effort. In order to suppress correlations in the x-dependence of the 
roughness, the length of each rectangular module is chosen at random in the range of 0.24nm (one unit cell) to lOnm. 
We then use the above Fibonacci-like procedure to assemble a scattering geometry [see figure [3|b)] of up to several 
jam in length. Finally, all modules are connected to two ideal half-infinite graphene waveguides. We average over 100 
realizations ^ of nanoribbons to eliminate non- generic features of particular ribbon configurations. 

In addition to the quantization steps due to transverse confinement [dashed black line in figure [sj^c)], a graphene 
nanoribbon with a perfect zigzag boundary of fixed width W features edge states with finite dispersion [figure |3|d)] 
since the coupling between the outermost carbon atoms is non-zero^. Consequently, the edge states of an ideal 
nanoribbon give rise to a peak in conductance G just below the Dirac point [shaded area in figure |3|c,d)]. In 
contrast to first-nearest-neighbor tight binding, and in line with the full ab-initio bandstructure and experimentP, 
our third nearest neighbor approach accounts for the breaking of electron-hole symmetry. Conductance is thus only 
approximately symmetric relative to E = 0. 

In the presence of edge disorder, G undergoes several pronounced changes [figure [sj^c)]: overall, G decreases with 
increasing distance in energy from the Dirac point relative to the ideal ribbon. The quantization steps due to 
the transverse confinement are strongly suppressed. Moreover, the edge disorder completely removes the sharp 
conductance peak attributed to edge states, as they are no longer conducting but become localized parallel to the 
ribboiP^, i.e. along the direction of transport. Consequently, corresponding signatures are difficult to observe in 
transport measurements of realistic samples. Scanning tunneling spectroscopy provides an alternative approach: 
peaks in the local density of states at energies slightly below the Dirac point have been recently, indeed, observed in 
STS experimentPl. 

In the limit where quantization steps due to the transverse confinement are strongly suppressed [figure [sjc)] pro- 
nounced broad dips in transmission (see arrows in figure |3| replace the original steps in conductance. This counter- 
intuitive reduction of transmission with increasing energy in the vicinity of steps can be qualitatively understood by 
considering Fermi's golden rule for the scattering of mode \nk) into mode \n^k^)^^^^ 

T{E) cc ^ I {nk\H\n^k^)f p^' (E). (18) 



Two trends contribute to this effect: firstly, strong fluctuations of the ribbon width broaden the DOS pn'{E) and 
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L[nm] VF[nm] 

FIG. 5. (a) Anderson-localized scattering state shown for a section (1100-1200nm) of a ribbon with total length L — 3/jm. 
(b) Conducting scattering state [as in figure |4]^a)] shown for the same section as in (a), (c) Longitudinal dependence 
of the conducting state [red, see (b)], and localized states [blue]. The latter are averaged over 100 Anderson-localized states 
for average ribbon width W = 20nm and different edge roughness amplitude AW (see insets), (d) Localization length I a as 
function of ribbon width W at energy E = 0.2eV for the five different values of edge roughness AW as in (c). 



smoother! the steps. Secondly, as the ribbon locally narrows, backscattering via scattering into evanescent modes is 
enhanced. This occurs preferentially for energies close to the opening of a new mode and results in a reduction in 
transmission causing the dips. 

It is instructive to compare the present results to calculations for edge- and bulk- disordered semiconductor nanowires 
featuring a parabolic dispersion relation in the long- wavelength (continuum) limit. While a reduction of quantization 
steps by disorder is observed for edge-disordered semiconductor ribbons^^ resembling the present results, the prominent 
transmission dips observed for graphene (arrows in figure 3) appear not to be present in such a system (compare with 
figure 2 in^^). However, dips have been found in other disordered semiconductor nanostructures that are associated 
with resonances supported by attractive (bulk) disorder potentialJ^. The distance in energy between these resonances 
and the quantization step (i.e. the subband minimum) corresponds to the binding energy of the (quasi) bound state^^. 
In the present case of graphene with rough edges, the enhancement of the local density of states near the Dirac point 
resulting from localized states at the edges is well knowiP^. Their statistical weight has been found to be much higher 
for graphene than for conventional nanostructureJ^. In the absence of a local attractive potential, these localized 
states could take on the role of resonances: for edge structures similar to the ones we investigate, the resonance 
energies Ei are statistically distributed in the range Ei G [—80, 0] meV^, and could, possibly, give rise to the observed 
broad dips. 



B. Localized scattering states at edges 

To gain a deeper understanding of the transport characteristics of edge-disordered graphene ribbons, we now analyze 
also individual scattering states. We find that states with energies where the conductance is only weakly perturbed 
by edge disorder [e.g. open triangle in figure |3jc)] feature a low amplitude at the edges [see figure ga)]. The overall 
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probability density of these scattering states remains concentrated near the center of the ribbon [see figure [4[a)] and 
is therefore only slightly affected by edge disorder. When only a single mode is open in the leads (at energies close to 
the Dirac point), modes located at K and K' in momentum space are not coupled in a zigzag graphene nanoribboiP. 
Only one cone contributes to transport in each direction [see figure [sj^d)]. This imbalance in the number of left- 
and right-moving channels on each cone is a special property of zigzag graphene nanoribbons^^, similar to the band 
structure of topological insulators. Backscattering is only possible in this energy window by inter-valley scattering at 
the rough edges. Since we observe a nearly perfectly conducting channel [figure |5|b),(c)], inter- valley scattering that 
requires momentum transfers of the order \K — K'\ is obviously suppressed at low energies^^. 

By contrast, scattering states at energies where transmission (and conductance) is considerably reduced [solid 
triangle in figure [sj^c) and figure [4]^b)] feature a strong enhancement of their wavefunction near corners of the edges 
originating from one sub-lattice only [see arrows in figure [4]^e, f)]. Projections onto the A and B sublattices [figure 
[4] (e, f)] show pronounced differences reflecting the violation of pseudo-spin conservation [equation We find 

enhancements of the A (B) sub-lattice scattering wave function at the upper (lower) edges of the ribbon, i.e., at 
those edges where the outermost carbon atom is of type A (B) [see zoom-ins in figure |4]^d)-(f)], in line with a strong 
enhancement of the local DOS near rough edge^^^^^. However, the pronounced differences in the wavefunction 
patterns near the center of the ribbon are not accounted for only by localized edge states since their decay length into 
the ribbon interior is much smaller than the ribbon width. We therefore attribute the dramatic drop in conductance 
to pronounced intra- valley and inter- valley backscattering at the edge corners, since the suppression of backscattering 
associated with the conservation of pseudo-spin [equation ^] no longer holds. 

As reported in earlier work on edge disorder in rough-edged graphene nanoribbons, transmission is strongly sup- 
pressed close to the Dirac point, leading to the formation of a transport gapl^^l^. Atomic-scale de fects on the edges of 
wide ribbons may lead to exponential (i.e., Anderson) localization due to destructive interferenc d"*^^ * ^^*^ -^. We use our 
modular approach to calculate scattering states on mesoscopic length scales [ribbon length L = 2/im, see figure [sja, 
b)]. By averaging over many realizations of edge disorder, we can thus explicitly probe for exponential localization 
and determine the localization length. Looking at the longitudinal dependence of the scattering state, 

nW 

mx)\'= / mx,y)\'dy, (19) 

we observe an exponential decay over up to 10 orders of magnitude [figure [sjc)]. Fitting to the functional form 

|'0(x)| oc exp(— x//a) we can numerically extract the localization length Ia- We find Ia to scale as Ia ^ aW/AW^ 
i.e., I A increases linearly with ribbon width and is inversely proportional to the disorder amplitude AW [figure |5|d)]. 
The localization length Ia is found to increase with increasing distance (in energy or in k) from the Dirac point (not 
shown), as suggested by the disorder-induced formation of a transport gapP^HI^. 

Superimposed on the exponential decay are oscillations on two shorter length scales: (i) a short beating period of 
A = 0.7nm due to interference between the K and cones^^ [A in figure |5|a)] and (ii) a much slower variation with 
the length scale A ^ 30nm [A in figure [sja)] which corresponds to the wavelength A = 27r//c associated with the linear 
dispersion relation E = v^hk^ i.e., the distance in k space from the K point. 

For comparison we also plot for the nearly perfectly conducting channel its wavefunction and its projection according 
to equation ([l9| [figure [5|b),(c)]. If the incoming scattering wave couples to the near-perfectly conducting channel, this 
contribution will be dominant after a certain ribbon length as all other contributions quickly die out. The oscillations 
due to K — K' interferences (i) are also present for this conducting state, though at reduced amplitude. While we 
observe Anderson localization for incoming scattering states at energies where more than one mode is open per cone, 
near-perfect conduction^^ appears to be confined to the topologically insulating part of the band structure. We 
expect these states to have a localization length that exceeds the dimension of our structure, if it is, at all, finite. 



C. Variations of edge roughness 



To investigate to what extent the above results depend on our particular choice of rectangular edge roughness 
we generalize our approach to include randomly jagged edges. We combine graphene segments featuring horizontal 
zigzag edges with segments featuring a boundary profile tilted by an angle /3 with respect to the horizontal zigzag 
direction [see insets in figure |6] . As outlined in section III , we calculate a set of modules (we use modules with length 
L G [20,40]nm) by direct inversion of a finite-sized Hamiltonian, and combine these modules to efficiently generate 
very long structures (total length > 1/im). 

Qualitatively, we find the same Anderson localization behaviour [see figure [oj^a)] as a function of ribbon width W 
and roughness amplitude AW (for fixed /3) as in the case with rectangular modules. However, unlike the case of 
free-particle dispersiorJ^, graphene nanostructures feature an interesting interplay between lattice orientation and 
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FIG. 6. (a) Localization length I a as a function of ribbon width W for different values of edge roughness amplitude AW (see 
inset). Each curve is averaged over 100 disorder configurations, featuring random edge directions, (b) Same as (a) for fixed 
AW = 2.0nm, for different edge roughness configurations: the parameter /3 labels the angle (in degrees) between straight- 
line segments and the horizontal zigzag direction of the graphene lattice (see top left inset), resulting in different roughness 
configurations (see bottom right inset). 



surface roughness. As this interplay is determined by the alignment angle /3 between the lattice orientation and the 
roughness, we can explicitly study its influence on transmission through the ribbon. We observe, indeed, that the 
value of the localization length strongly depends on the shape of the boundary with respect to the discrete lattice: 
edges consisting of randomly concatenated zigzag-edges only (i.e., with /3 = 60°) show substantially longer localization 
lengths than edges formed by an even mixture of zigzag and armchair edges [i.e., with /3 = 75°]. The dependence 
of the localization length on /3 can be understood in terms of the length of undisturbed zigzag (or armchair) edges: 
cutting close to a symmetry plane of the lattice (i.e. 60° or 30°) results in comparatively longer segments of zigzag (or 
armchair) edges. By contrast, a cut at 75° yields an irregular sequence of very short segments of armchair and zigzag 
boundaries and thus strongly breaks the translation symmetry of a clean zigzag (or armchair) edge. As an aside we 
note that the data presented in the previous subsection includes a variation in ribbon direction, i.e., the ribbons are 
not perfectly straight, notably in figure [sja, b). This amounts to an effective increase of edge roughness. As a result, 
the localization length is further decreased in that case [compare figure [5|d) with figure |6|b) for (3 = 90°]. 

The observation of the relative change in resistance as a function of /3, i.e., of the angle between the graphene lattice 
and the atomic-scale edge, might have implications for experiments. Measuring the atomic-scale roughness is difficult 
requiring an STM setup. Measuring localization length for different ribbon widths might provide an alternative probe 
for the atomic-scale edge roughness. Conversely, our results could be tested by comparing transport measurements 
for nanoribbons fabricated with different methods (i.e. etching, growth on Si-C substrate^^, unzipping of graphene 
nanotubes^^) resulting in (known) different edge characteristics. 



D. Fourier analysis of channel states 

We explore now the interplay between short-range defects in real space and the absence (presence) oi K — K' 
inter-cone scattering in /c-space. For this purpose we analyze the Fourier transforms of the asymptotic scattering state 
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transmitted 




K' K K' K 



FIG. 7. Two-dimensional Fourier transform |'0(k)|^ [top row], and longitudinal dependence [bottom row, equation 

(24)] of the incoming [(a,d), left column], reflected [(b,e), center column] and transmitted [(c,f), right column] part of the 
scattering state in the waveguides. Due to the finite size of the numerically evaluated scattering state, the Fourier transform 
features a grid of thin horizontal and vertical lines. The insets show an enlarged view of the point [a dashed white circle is 
inserted as guide to the eye]. The first Brillouin zone of the reciprocal lattice is shown as white hexagon. 



in the semi- infinite entrance (exit) waveguides. The Fourier transform is calculated as 



V^(k) = T[iIj{t)] = / dr V^(r) e 
J A 



i kr 



(20) 



where we extend the integral over a finite area A in the asymptotic region of the waveguide, i.e., far away from the 
scattering region. Three different classes of asymptotic scattering states need to be considered: (i) the incoming Bloch 
states propagating in x-direction with wavenumber /c^, 



Xniy), 



(21) 



where Xn{y) represents the transverse eigenfunction of mode n of the semi-infinite nanoribbon, (ii) the waves trans- 
mitted through the disordered region with transmission amplitude tmn^ and (iii) the reflected waves with reflection 
amplitude rnm- The corresponding Fourier components are given by 



(22) 



(23) 



For a hexagonal lattice, the real and the reciprocal lattice are rotated by 90 degrees with respect to each otheil^. 
The first Brillouin zone for the ideal zigzag ribbon is thus given by a hexagon resting on a side rather than on a tip 
[see white hexagon in figure [t]^ a)]. To better visualize the enhancement of '^(k) near the K or K' points, we integrate 
i^Ck) over the transverse direction. 



I 



dky\m)\' 



(24) 
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In perfect zigzag ribbons, incoming modes feature non- vanishing amplitudes either near the K or the K' point, i.e., 
there is no couphng (scattering) between K and K' . For the incoming Bloch state the projected Fourier transform 
equation ( [24| features peaks at the kx values corresponding to K' [see figure [rj^d)]. The close-up of the peak in V^(k) 
near the .S^^point [inset of figure [7|^a)] is structureless for the incoming Bloch wave with fixed transverse quantum 
number n. The horizontal and vertical lines are finite-size effects of the Fourier-transformed sample. Likewise, the 
origin of the additional bright spots inside the Brillouin zone is zone folding (the Brillouin zone of the ribbon is smaller 
than the graphene Brillouin zone). The interesting physics, on the other hand, is contained in finite amplitudes at 
both K and K' points of the scattered wave [see figure [7]^b,c)] which are induced by i^T—i^' scattering at rough edges. 
The relative strength of the integrated K and K' peaks [see figure [7|^e,f)] is a direct measure for the amount of 
inter-valley scattering. Furthermore, we observe a pronounced fine structure near the K' and K points: enhancement 
along a half-circle forms around the Dirac points [see inset in figure [7|^b,c)]. The surface of section of the double-cone 
band structure of constant energy is approximately a circle, the diameter of which is proportional to the energy. 
In the refiected (transmitted) part of the wavefunction, we only see the left half (right half) of this circle being 
populated, corresponding to negative (positive) group velocities. Enhancement along the full semicircle is due to 
inter- mode scattering n ^ m between transverse modes within the same valley (intra- valley scattering). We can thus 
conclude that pronounced intra-valley scattering at the rough edges distributes the refiected (or transmitted) wave 
almost uniformly over the energetically accessible half-circle of the band structure compatible with their propagation 
direction. The Fourier transform thus allows us to assess the amount of both inter- valley K — K' scattering (by the 
relative amplitude around the K and K' points in the reciprocal lattice) and inter-mode scattering by the angular 
distribution on the half-circle of a single cone for the incoming and reflected (or transmitted) states. 



V. CONCLUSIONS AND OUTLOOK 



We have presented a novel numerical approach to efliciently calculate the Green's function for extended nanoribbons. 
Key is the build-up of the ribbon by a random assembly of modules. Connecting these modules by way of a Dyson 
equation allows us to calculate the transport properties of long graphene ribbons. We find the conductance to be 
suppressed by rough edges relative to that of the perfect ribbon. For low energies, we observe near-perfectly conducting 
channels due to the band structure of zigzag graphene nanoribbons. Quantization steps are washed out and, in part, 
replaced by dips due to scattering into evanescent modes'^, in contrast to edge-disordered semiconductor nanoribbons 
with free-particle dispersion^^. An analysis of individual scattering states in both real space and Fourier space reveals 
pronounced A — B sublattice asymmetries and K — K' scattering. We determined specific signatures of inter- and 
intra-valley scattering by Fourier transform spectroscopy of scattering states. We also identified Anderson localized 
states for different disorder configurations, extending over several micrometers with an exponential decay spanning 10 
orders of magnitude. The corresponding localization length was calculated as a function of both the magnitude of edge 
roughness, and its alignment with the graphene lattice. We find that the latter plays a significant role in determining 
the localization length hinting at the importance of correctly modeling microscopic details of edge disorder beyond 
its amplitude and correlation length. 

We conclude by pointing to possible future applications. While early transport measurements were strongly affected 
by substrate interactions resulting in puddles of electron and hole conductivity due to bulk disorder, recent advances in 
the manufacturing of much cleaner graphene nanostructures by growing on Si-C substrates^^, unzipping nanotubes to 
arrive at smooth-edged ribbonJ^ as well as suspended graphene^ have shifted the focus to edge disorder investigated 
in the present work. Indeed, the measurement of size quantization plateaus has been surprisingly elusive in graphene 
nanoribbons^SMH^ in qualitative agreement with our present findings. Only recently, by prolonged annealing and 
suspending graphene nanoconstrictions, first signatures of size quantization could be found^*". Our findings regarding 
the sublattice sensitivity of the wavefunction (figure 3) could be tested by STM scans of bound states in graphene 
nano-islands^^: in these measurements strong enhancements of wavefunction amplitudes were found on one sublattice, 
resulting in trigonal patterns close to edges that are quite similar to our numerical findings. Our simulations predict 
similar STM patterns for scattering states in extended nanostructures. In particular, such measurements could 
elucidate the precise nature of the scattering mechanisms encountered at edges prepared wth different techniques. 
Indeed, we expect defects affecting both sublattices, as recently investigated^^, to exhibit different signatures than 
e.g. single vacancies. Finally, magnetic field effects allow for an additional external parameter more easily tunable in 
experiments than lattice geometries. States associated with different K points react differently to magnetic fields. This 
dependence might help to disentangle contributions from different i^-points in scattering accessible by our Fourier 
analysis. We note that our algorithm can be easily adapted to accomodate magnetic fields while retaining its favorable 
scaling properties. Investigations in this direction are currently under way. 
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